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Abstract. - We investigate global persistence properties for the non-equilibrium critical 
dynamics of the randomly diluted Ising model. The disorder averaged persistence probability 
Pc{t) of the global magnetization is found to decay algebraically with an exponent 9^ that we 
compute analytically in a dimensional expansion in d = 4 — e. Corrections to Markov process 
are found to occur already at one loop order and 6c is thus a novel exponent characterizing 
this disordered critical point. Our result is thoroughly compared with Monte Carlo simulations 
in d = 3, which also include a measurement of the initial slip exponent. Taking carefully 
into account corrections to scaling, 9c is found to be a universal exponent, independent of the 
dilution factor p along the critical line at Tc{p), and in good agreement with our one loop 
calculation. 



Persistence properties, and related first-passage-time problems, have been intensively in- 
vestigated both theoretically and experimentally these last few years [1]. Originally introduced 
in the context of non-equilibrium coarsening dynamics of ferromagnets at zero temperature [2] , 
the persistence exponent describes the long time decay P(t) ~ of the probability P{t) 
that a stochastic variable, here the local order parameter, does not cross a threshold value 
during a time interval t. Interestingly, this exponent 9 appeared to be highly non trivial, even 
for simple pure systems such as diffusion [3], and therefore it has motivated numerous inves- 
tigations in a various pure models. However, except for few models in one dimension [4,5], 
very little is known about this exponent for disordered systems. 

Although for a ferromagnet at finite T, one expects that the persistence P{t) associated 
to the local magnetization has an exponential tail, persistence properties received a special 
attention in the context of critical dynamics of pure ferromagnets at Tc- Indeed, it has been 
proposed [6] that the probability Pc{t) that the global magnetization M has not changed 
sign in the time interval t following a quench from a random initial configuration, decays 
algebraically at large time Pc{t) ^ t"^". In this context, analytical progress is made possible 
thanks to the property that, in the thermodynamic limit, the global order parameter remains 
Gaussian at all finite times t. Indeed, for a d-dimensional system of linear size L, M{t) is 
the sum of random variables which are correlated only over a finite correlation length 
£,{t). Thus, in the thermodynamic limit L/£^{t) 3> 1, the Central Limit Theorem (CLT) 
asserts that M{t) is a Gaussian process, for which powerful tools have been developed to 
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compute the persistence properties [1,7]. Remarkably, under the additional assumption that 
M is a Markovian process, 9c can be related to the other critical exponents via the scaling 
relation 0c = /^=(A — d+l — ri/2)z^^ , with z and A the dynamical and autocorrelation [8, 9] 
exponent respectively. However, as argued in Ref. [6], M is in general non Markovian and 
thus Oc is a new exponent associated to critical dynamics. For the non conserved critical 
dynamics of pure O(A^) model, corrections to this scaling relation were indeed found to occur 
at two loops order [10], in rather good agreement with numerical simulations in dimensions 
d = 2, 3 [11]. The question whether such a wnwersa^ persistence exponent can also be defined, 
and computed, in the presence of quenched disorder remains an open question. 

On the other hand, there is currently a wide interest in the slow relaxational dynamics 
following a quench at a critical point [12]. Although simpler to study than glasses, they 
display interesting properties such as two time dynamics (aging) or violation of the Fluctuation 
Dissipation Theorem, as found in more complex glassy systems [13]. Recently, some progress 
has been achieved in the characterization of the effects of the disorder on these properties [14- 
16]. In this Letter, taking advantage of these recent studies, we study persistence properties, 
both analytically using RG and numerically via Monte Carlo simulations, of the randomly 
diluted Ising model: 

H^^p^PjS.Sj (1) 

where Si are Ising spins on a d-dimensional hypercubic lattice and pi are quenched random 
variables such that pi ~ 1 with probability p and otherwise. For the experimentally relevant 
case of dimension d = 3 [17], for which the specific heat exponent of the pure model is positive, 
the disorder is expected, according to Harris criterion [18], to modify the universality class of 
the transition. For 1 — p <C 1, the large scale properties of 1^ at criticality are then described 
by the following 0(1) model with a random mass term, the so-called Random Ising Model 
(RIM): 



l(V^)2 + i[r„ + ^(x)]^2 + |^4 



(2) 



where i/j = ^p{x) and i'{x) is a Gaussian random variable '4i{x)'4i{x') = Ao(5'^(a; — x') and tq, 
the bare mass, is adjusted so that the renormalized one is zero. 

We first study the relaxational dynamics of the Randomly diluted Ising Model 10) in 
dimension d = A — e described by a Langevin equation: 

where {C{x,t)) = and {(^{x,t)(^{x' ,t')) — 2rjT5{x — x')5{t — t') is the thermal noise and 77 the 
friction coefficient, set to 1 in the following. At initial time U — 0, the system is in a random 
configuration with zero mean magnetization and short range correlations [Lp{x^t — 0)ip(y,t = 
0)]i = T^^S'^{x — y), where Tq"^ is irrelevant in the RG sense [8] (it will thus be set to zero in 
the following). Defining the global magnetization M{t) as 

Mit) = ^^As,(t) = / ^ix,t) (4) 

-'^occ ^ ^ Jx 

where iVocc is the total number of occupied sites, we are interested in the disorder averaged 
probability Pc(t) that the magnetization has not changed sign in the time interval t following 
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the quench. In that purpose, we introduce the correlation Cu^ — {M{t)M{tw)) and the hnear 
response Ti-tt^ to a small external field f{tw), TZu^ = 5{M{t)) / 5f{tw) where — and (...) denote 
respectively averages over the disorder and thermal fluctuations. At one loop order, Ctt„ and 
TZtt^ take the scaling forms [15], Cu^ — A^t^ ^(.(t/tw) and 7^ft„ = Aj.t ^<i>r(V^t") where 
Ac, Ar are non universal constants and, ^c.r{x) - not shown here for clarity - are universal 
scaling functions such that <I>c,r(l) = 1- In a previous publication [16], we have shown, using 
the Exact Renormalization Group for the dynamical effective action [19], that they are given, 
at this order, by the solution of the following equations 

dtUtt^ + fi{t)ntt^ + dtiY.u,nt,t^ = , ^(i) = - ^ dh^u, (5) 
c,,^ = 2T / " dhUtt^TZt^t^ + f dh r dhn,t^Dt,t,nt^,^ (6) 

Jti Jti Jti 

where the self energy Ef^ta and the noise-disorder kernel I?tit2 have been computed in [16] : 

^t,t,^-l^{l{ti~t2)f , Alt. = yy^(7(ii-<2)- 7(^1+^2)) (7) 

where 7(2:) = {x + A(^^)~^, Aq being the UV cutoff [20]. These equations Eq. lOlHI), suggest 
to describe the time evolution of the magnetization by an effective Gaussian process M{t): 

dtM{t) + ^l{t)M{t) = - / dh^mMih) + Cit) (8) 
Jo 

where (^{t) is an effective disorder induced Gaussian noise with zero mean and correlations 
(C(OC(^'))off — 2TriS{t — t') + Dtt'- The idea is then to compute Pc{t) as the persistence 
probability of the process M{t). The rhs of Eq. JSJ clearly indicates that this process is non- 
Markovian. However, as Stf as well as Du' are of order ©(y^), one can use the perturbative 
computation of 9c around a Markov process initially developed in Ref. [7] and further studied in 
the context of critical dynamics in Ref. [10] . In that purpose [1] , let us introduce the normalized 



Gaussian process m{t) — M {t)/ \J {M^{t))cs. Let T = Ini, then m(T) is a stationary Gaussian 
process, its persistence properties are then obtained from its autocorrelation function: 

{m{T)m{T^)),B = e-'^(^-^™U(e^-^™ ) (9) 

Under this form (Q, one can use the first order perturbation theory result of Ref. [10] to 
obtain the one loop estimate: 




;z=W— — + 0(e) = 0.06968... in d = 3 (10) 



where is the value corresponding to a Markov process. The second term in Eq. (|10|l is the first 
correction due to the non-Markovian nature of the dynamics. Interestingly, this correction is 
entirely determined by the non-trivial structure of the scaling function A{x) , which is directly 
related to ^c{x). Notice also that, at variance with the pure 0{N) models [6, 10], these 
corrections in presence of quenched static disorder already appear at one loop order. 
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We now turn to the results from our Monte Carlo simulations of the relaxational dynamics 
of the randomly diluted Ising model ^ in dimension d — 3, which were done on cubic 
lattices with periodic boundary conditions. For p = 1, the system is pure and shows ferro- 
magnetic order at T < Tc{p — 1). The critical temperature Tc{p) [21,22] decreases with p 
and becomes at the percolation threshold p = Pc (for a c? = 3 cubic lattice, pc ~ 0.341). 
For p > Pc, the system is initially prepared in a random initial configuration with zero mean 
magnetization Mq = 0. Up and down spins are randomly distributed on the occupied sites, 
mimicking a high-temperature disordered configuration before the quench. At each time step, 
one site is randomly chosen and the move Si — > —Si is accepted or rejected according to 
Metropolis rule. One time unit corresponds to such time steps. The exponent dc is mea- 
sured numerically for cubic lattices of linear size L = 8, 16, 32 and 64. After a quench to Tc 
from the initial random configuration each system evolves until the global magnetization first 
change sign. Pc{t) is then measured as the fraction of surviving systems at each time t, over 
a number of samples which varies from 2 x 10^ for L = 8 to 2 x 10'* for L = 64. Although 
the Renormalization Group analysis stays that the critical exponents are independent of p (at 
least for 1 — p <gC 1) the question of universality in disordered systems, and in particular in the 
present problem [21], is a longstanding issue. Therefore we will compute Pc{t) for different 
values of p{> pc) — 0.499, 0.6, 0.65 and 0.8 along the critical hne at Tc{p). 

In fig.^we present the results of Pc{t) for p — 0.8 and for different lattice sizes. According 
to standard finite-size scaling [6], one expects the scaling form Pc{t) = t^^" f{t/L^), where z 
is the dynamical exponent. Keeping the rather well established value of z = 2.62(7) [16,23] 
fixed, 9c is varied to obtain the best data collapse. The final scaled plot is shown in the inset 
of fig. ni We have checked (not shown) that the same analysis can be done for all other values 
of p studied here. However, the best data collapse is then obtained for a value of 9c which is p- 
dependent. This is depicted in fig. El (left) where one clearly observes a decrease in slope while 
occupation probability p is reduced, which could naively suggest a non universal value for the 
exponent 6c{p)- However, it is now well known that the present model is strongly affected by 
corrections to scaling [16,23], characterized by a universal exponent 6, i.e. independent of p, 
which has been estimated in previous works [16,23] to be & = 0.23(2). Therefore following 
these previous studies we suggest that the persistence can be written as 



Pcit)^t-'>^Mt) , fp{t)^A{p)il + B{p)t-') 



(11) 
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Fig. 2 - Left: Persistence probability Pc{t) for different p = 0.5, 0.6, 0.65 and 0.8 with fixed L — 64. 
The system is initially prepared with A/q ~ 0. Inset: Universality of Pc{t) for different p. The 
function fp{t) is defined in the Eq. lllll . Right: Similar to fig. |21(left), but prepared with a non-zero 
initial magnetization AIo = 0.001. The linear system size in this case is L = 100. 



which A{p) and B{p) arc fitting parameters. As shown in the inset of fig. |21 (left), the good 
collapse of the quantity P{t)/ fp{t) for different values of p suggest that 9^ is a universal 
exponent. A{p) and B{p) are found to be monotonous decreasing function of p and such that 
B{p — 0.8) ~ 0. The value extracted from this collapse is 

6lc = 0.35 ±0.01 (12) 

We have also computed the persistence probability for systems quenched from random 
configuration with a small initial magnetization Mq = 0.001. The number of up A^up and 
down A'down spins are thus : A'up = (1 + Mo)/2Nocc and iVdown = Aocc — N^p. First we 
randomly distribute the A'^^p up spins in the occupied sites of the lattice and then fill up the 
rest with down spins. As noticed previously [11] this protocol allows to reduce the statistical 
noise and thus to study larger system sizes (this however renders the finite size scaling analysis 
more subtle [11]). In fig. |21 (right), we plot the persistence for system size L — 100 ( the data 
have been averaged over 2 x 10** ensembles) for p — 0.499, 0.6, 0.65 and 0.8. The study of 
larger system size allows to reduce the corrections to scaling. Indeed, although in the short 
time scales, the straight lines have slightly different slopes which depends upon p, at later 
times the slopes varies from 0.36(1) for p — 0.8 to 0.35(1) for p — 0.499 : this confirms the 
p-independent value of dc obtained previously Eq. 1)12(1 . 

In order to compare this numerical value ((12|l with our one loop calculation pO(l one 
needs an estimate for ^. Such an estimate is needed not only for the sake of this comparison 
but also to characterize quantitatively non Markovian effects. The argument mentioned in 
the introduction, relying on the CLT, which stays that the global magnetization is, in the 
thermodynamic limit, a Gaussian variable is also valid in the presence of disorder, and we have 
checked it numerically. Therefore, a finite difference A 1(10(1 is the signature of a Non Markovian 
process. A numerical estimate of /i can be done in two different ways [24]. First, using the 
numerical estimates for t] = 0.0374(45) from Ref. [22], z = 2.62(7) from Ref. [23], and X/z = 
1.05(3) from Ref. [16] one obtains /inum = 0.28(4), which gives a first estimate Anum = 0.07(5). 
Because of the relatively big error bars on A and z, we propose, alternatively, to express in 
terms of the initial slip exponent , 9' [8], using X — d— z9' , a,s ^ ~ — 0' + (1 — ry/2)/z. Although 
9' has been estimated numerically for pure Ising systems [25,26], there are no available data 
for the disordered case : we thus now present a numerical computation of 9' . 
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Fig. 3 - Rescaled global magnetization M{t) / gp{t) as a function of t for p = 0.49, 0.6, 0.65 and 0.8 
in the log-log scale. The linear system size is L = 64, the initial magnetization is Mq = 0.01 and the 
measured exponent 6' — 0.1. Inset: Global magnetization M{t) for p = 0.8 and L — 16, 32 and 64. 



In that purpose we study the time evolution of M{t) when the system is quenched from 
an initial configuration with short range correlations but a finite, however small, magneti- 
zation Mq. The initial stage of the dynamics is characterized by an increase of the global 
magnetization, described by a universal power law [8] 

M{t) - Mot"' (13) 

At larger times, t ^ to, M{t) decreases to zero as M{t) - In our simulations, the 

system is initially prepared in a random initial configuration with a small but finite initial 
magnetization Mo = 0.01 and quenched down to Tc at t = 0. In all subsequent times we 
measure M{t) Q for linear system sizes L — 8, 16, 32 and 64. Finally data are averaged over 
8x 10^ samples for L = 8 to 10^ samples for L — 64. In the inset of fig. El we show a plot of M{t) 
for p = 0.8 and different system sizes. One sees clearly that M{t) is increasing until a time to, 
which is an increasing function of the system sizes considered here [27], and then decreases 
to zero (although the aforementioned scaling for t ^ to is not clearly seen here). Here also, 
we compute M{t) for different values of p = 0.499, 0.6, 0.65 and 0.8 and observe corrections 
to scaling which we take into account as M{t) = t" gp{t) with gp{t) = A'{p)[l + B'{p)t~^], 
A'{p), B'{p) being fitting parameters. As shown in fig.|31 one obtains a reasonably good data 
collapse of M{t)/gp{t) vs. t for the different values oi p. Here also, one finds that p — 0.8 
is almost unaffected by these corrections to scaling, i.e. i3'(0.8) ~ 0. After a microscopic 
time scale, which increases as p is lowered, one observes a universal power law increase (|13|l . 
from which we get the estimate 0' ~ 0.10(2). This value is in good agreement with the 
two-loops estimate ^2 loops = 0.0868 [28] and agrees also quite well with the scaling relation 
A = d — z6' [16]. This gives //num = 0.27(3) and our final numerical estimate 

Anum = 0.08 ± 0.04 (14) 

which is in good agreement with our previous one loop estimate (|10|) . We also notice that 
these deviations from a Markov process are slightly larger than for the pure case [11]. 

In summary we have shown that the global persistence exponent 9c carries the signature of 
Non Markovian effects, and is thus a new critical exponent characterizing this random critical 
point, that we have computed within a one loop approximation ()10|l. Our detailed numerical 
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analysis, which rehes upon a computation of 6^ and the initial slip exponent 9' , supports the 
existence of these Non-Markovian violations (|14|l and is in good quantitative agreement with 
our perturbative approach. In addition, Oc is found to be universal along the critical line. In 
view of recent progress in the study of aging properties in finite dimensional glassy phases [29] , 
it would be very interesting to extend the approach presented here to these situations. 
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